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Abstract —With the increasing adoption of plug-in electric vehi¬ 
cles (PEVs), it is critical to develop efficient charging coordination 
mechanisms that minimize the cost and impact of PEV integra¬ 
tion to the povrer grid. In this paper, we consider the optimal PEV 
charging scheduling, where the non-causal information about 
future PEV arrivals is not known in advance, but its statistical 
information can be estimated. This leads to an “online” charging 
scheduling problem that is naturally formulated as a finite- 
horizon dynamic programming with continuous state space and 
action space. To avoid the prohibitively high complexity of solving 
such a dynamic programming problem, we provide a Model 
Predictive Control (MPC) based algorithm with computational 
complexity 0{T^), where T is the total number of time stages. 
We rigorously analyze the performance gap between the near- 
optimal solution of the MPC-based approach and the optimal 
solution for any distributions of exogenous random variables. 
Furthermore, our rigorous analysis shows that when the random 
process describing the arrival of charging demands is first-order 
periodic, the complexity of proposed algorithm can be reduced 
to 0(1), which is Independent of T. Extensive simulations show 
that the proposed online algorithm performs very closely to the 
optimal online algorithm. The performance gap is smaller than 
0.4% in most cases. 


I. Introduction 
A. Background and Contributions 

The massive deployment of PEVs imposes great challenges 
to smart power grid, such as voltage deviation, increased power 
losses, and higher peak load demands. It is critical to design 
PEV charging mechanisms that minimize the cost and impact 
of PEV integration. Previously, PEV charging coordination has 
been extensively studied to minimize power loss, minimize 
load variance, or minimize charging cost, etc 0-0- Ideally, 
the load demand can be flattened as much as possible if 
the information about future charging demand is known non- 
causally when calculating the charging schedule. However, in 
practice, a PEV charging station only knows the load demand 
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of the PEVs that have arrived, but not that of the PEVs 
coming in the future. Eortunately, the statistical information 
of the future charging demands can often be acquired through 
historic data, which benefits the control of the PEV charging 
scheduling in practical scenarios. 

In this paper, we consider the optimal PEV charging 
scheduling, assuming that the future charging demand is not 
known a priori, but its statistical information can be estimated. 
In particular, we define the cost of PEV charging as a general 
strictly convex increasing function of the instantaneous load 
demand. Minimizing such a cost leads to a flattened load de¬ 
mand, which is highly desirable for many reasons 0-0- The 
online PEV charging scheduling problem is formulated as a 
finite-horizon dynamic programming problem with continuous 
state space and action space. To avoid the prohibitively high 
complexity of solving such a dynamic programming problem, 
we provide a Model Predictive Control (MPC) approach to 
obtain a near-optimal solution. Instead of adopting the generic 
convex optimization algorithms to solve the problem, we 
propose an algorithm with computational complexity 0{T^) 
by exploring the load flattening feature of the solution, where 
T is the total number of time stages. We rigorously analyze 
the performance gap between the near-optimal solution of the 
MPC-based approach and the optimal solution, and the result 
applies to any distributions of exogenous random variables. 
Specially, the performance gap is evaluated by the Value 
of the Stochastic Solution (VSS), which represents the gap 
between the solution of the approximate approach and that 
of dynamic programming problem Q-| 101. Eurthermore, our 
analysis shows that when the random process describing 
the arrival of charging demands is first-order periodic, the 
complexity of proposed algorithm can be reduced to 0(1), 
which is independent of T. Extensive simulations show that 
the proposed algorithm performs very closely to the optimal 
solution. The performance gap is smaller than 0.4% in most 
cases. As such, the proposed algorithm is very appealing for 
practical implementation due to its scalable computational 
complexity and close to optimal performance. 

The rest of the paper is organized as follows. A review of the 
related work on the PEV charging scheduling with uncertain 
load demand is presented in Section |I-B| We introduce the 
problem formulations of both offline and online PEV charging 
problem in Section |I^ In Section |I^ we propose a MPC based 
online algorithm and analyze its performance gap. The 0(1)- 



complexity algorithm is given when the arrival process is first- 
order periodic in Section IV Simulation results are presented 
in Section Finally, the paper is concluded in Section VI 


B. Related Work 

The works on the PEV charging scheduling with uncertain 
PEV load demand include both simulation-based evaluations 
lE)^ HI and theoretical performance guarantees 
Meanwhile, MPC is one of most commonly approaches for 
which has been widely adopted in recent studies |Tg-|Tg. 
© leverages the MPC based method to design a dynamic 
charging and driving cost control scheme. Both fT?) and |14| 
apply MPC algorithms to minimize the load variation. |15| 
proposes a plug and play MPC approach to minimize the 
voltage fluctuations by assuming that the load demand is time- 
periodic. Compared to Gg-in), in this paper we analyze 
the performance gap between the solution of MPC approach 
and the optimal solution regardless of the distribution of the 
load demand. Besides, we provide a more scalable algorithm 
with 0(1)— complexity as well as the optimality analysis 
for the case when the load demand is first-order periodic. 
Additionally, the objective functions in d-Gl) are quadratic 
forms of load demand. Whereas in this paper, the objective 
function is a general strictly convex increasing function which 
reflects both the charging cost and the load variance. 

As to the amount of information needed, the EV charging 
scheduling algorithms in and GZ) require the probability 
distribution of the random PEV arrival process. In contrast, 
the proposed algorithm in this paper only requires the first- 
order moment, i.e., the expected values of the random demand 
patterns. In practice, it is a lot easier to obtain the expected 
values of random process than to obtain the probability dis¬ 
tribution of a random process. Convex cost functions are also 
considered in ID and GD- Both of them devise the online 
algorithms for battery charging control problems, where there 
is no charging deadline for the battery. The PEV charging 
scheduling problem in this paper differs from stationary battery 
charging in that each PEV has a demand to be satisfied before 
a certain deadline. 


II. Problem Eormulation 

We consider the PEV charging scheduling problem, where 
PEVs arrive at the charging station at random instants with 
random charging demands that must be fulfilled before a 
random departure time. 

A. Optimal Offline PEV Charging Problem 

Eor the ease of understanding, we first introduce an ideal 
case, where all the non-causal information of based load and 
the PEVs, including the arrival times, departure times and 
charging demands are known to the charging station before 
the system time. The entire system time is divided into T 
equal-length time slots. Let Af denote the set of PEVs that 
arrive during the system time. Notice that for a given time 
slot number T, Af is itself a random set due to the random 
arrival of PEVs. We denote by 2{t) the set of PEVs that 


(s') 

are in the charging station during slot t. Denote by t) ’ and 
the arrival time slot and the departure time slot of PEV 
i, respectively, di denotes the charging demand that PEV i 
requires. The charging station needs to decide, the charging 
rate Xit,yi G Hit). To satisfy the demand di, xn must satisfy 

j(s) = di- Let St be the total charging rate of time slot 
t, i.e., 

st= ^ Xit,yt = 1,2,-■■ ,T, (1) 


which is also called charging load at time t. The total load 
consists of both the charging load and the inelastic base load in 
the same location. The base load, denoted by It, represents the 
load of other electricity consumptions at time t except for PEV 
charging. Then, the total load at time t is given by 
It- Suppose that the charging cost at time f is a strictly convex 
increasing function of the total load, denoted by /(s* -I- It). 
The convexity and increasing property of /(s* -I- /*) reflects 
the fact that each unit of additional power demand becomes 
more expensive to obtain and make available to the consumer. 
Eor example, in the wholesale market, the instantaneous cost 
can be modeled as an increasing quadratic function of the 
instant load Q-® . On the other hand, the convexity of f{st + 
It) also capmreOhe intent of reducing the load fluctuation 
over time |7). Then the total cost over time T is computed 
as fi^t + It). In the ideal case, assume that lt,t[^\t[^\ 

and di for all f = 1, • • • ,T,i G AA sre known non-causally at 
the beginning of the system time. Then, the charging station 
can solve •G) and obtain the optimal charging rate, denoted by 
for all time t and the optimal total cost, denoted by T*!. 
Such a solution is referred to as an “optimal offline solution”. 


tlii =min 


s.t. 


t=i \iex(t) 


yZ Xit = di,Vi G Af, 


*« > 0, Vt = , • ■ • , ^,yi GAf. 


(2a) 


(2b) 

(2c) 


In particular, the optimal total charging rate, denoted by 
S(, is defined as sj" = J^ieiiti^it- Note that there are in 
total 0(T\X{t)\) variables in Q, where \X{t)\ denotes the 
cardinality of the set X{t). This number can be quite large 
when the number of cars present at each time slot, |T(f)|, is 
large. Next, we propose an equivalent transformation of Q 
that drastically reduces the number of variables. In particular, 
the following Theorem [T] shows that as long as we find the 
optimal St Vf, the optimal x*f Vi, t can be obtained by earliest 
deadline first (EDF) scheduling. 


Theorem 1: If a set of stS satisfy the following inequality 
for all n = 1, • • • , T 

Y, Y di<^st<^ yz di, (3) 

then there exists at least a set of Xi^s that is feasible to (0. 
One such set of xu’s can be obtained by EDE scheduling, 
which charges the PEV i G X{t) with the earliest deadline 





at a rate St at each time t. Moreover, when St = s^, the set 
of Xit’s obtained by EDF scheduling are the optimal solution, 
x*t, to 0. 

Proof: Please see the detailed proof in Appendix IM To 
see Theorem note that 0 implies that the total energy 
charged by any time slot n is no less than the total charging 
demand that must be satisfied by time n and no more then 
the total charging demand of PEVs which have arrived up 
to time n. On the other hand, by EDE scheduling, PEVs 
with earlier deadlines must be fully charged before those 
with later deadlines can be charged. Thus, 0 guarantees 
the fulfillment of the charging demands of each individual 
PEV. With Theorem we can transform 0 to the following 
equivalent problem with T variables. 


T 

VPi^min + (4a) 

•St 

t = l 


n n 


s.t. 

E- 

>E E di,Vn = l,-- 

• ,T. 

(4b) 


t=i 





n 

E- 

n 

<E E di,Vn=l,-- 

■ ,T. 

(4c) 


t=i 





y^t=i 
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The optimal solution to 0 has an interesting feature: it 
does not change with the cost function f{st + h), as long as 
/ is strictly convex. Moreover, also minimizes the variance 
of total load subjecting to ( |4bl l and ( |4cl l, where the variance 

I'T’ 

of total load is defined as X]t=i('®t + h — 
m This is proved in Theorem]^ 

Theorem 2: The optimal solution Sj to 0 does not change 
with the cost function /(.), as long as /(.) is strictly convex. 
Moreover, is a load flattening solution that minimizes the 
variance of total load. 

Proof: Please see the detailed proof in Appendix 

Remark 1: In practice, a capacity constraint on s* + f is 
present for each t due to the hardware limitations and security 
concerns. The constraint is omitted in our formulation for 
the following reason. Theorem indicates that the optimal 
solution sj' to 0 minimizes the variance of total load. That is, 
any other scheduling solution would have a higher peak load, 
and therefore is more likely to violate the capacity constraint. 
In this sense, the optimal solution sj" to 0 is “capacity 
optimal” in the sense that if the optimal solution to Problem 
12 (or equivalently Problem 0 violates the capacity constraint, 
then there does not exist any other scheduling solutions that 
satisfy the capacity constraint. 


B. Online PEV Charging problem 

For the online PEV charging problem, the charging sched¬ 
ule only depends on the statistic information of future load 
demand, the current based load and the remaining charging 
demands and deadlines of the PEVs that have arrived so far. 
In contrast to the ideal case in the last subsection, in practice 
the charging station only knows the remaining charging de¬ 
mands and departure deadlines of the PEVs that have already 
arrived, as well as the current base load. The online charging 
scheduling algorithm computes the charging rate Sk at each 
time slot k based on the known causal information and the 


statistical information of the unknown future demands. The 
charging rate Sk, once determined, cannot be changed in the 
future. Specifically, the remaining charging demand of PEV 
i at time k is given by = di — Xa- Note that, 

d^ = di for all PEVs that have not yet arrived by time k — 1. 
A close look at 0 suggests that the charging schedule st 
only depends on the total charging demand that needs to be 
finished before a certain time, but not the demand due to 
individual PEVs. Thus, for notational simplicity, we define 
df = = k, - ■ ■ ,T, as the total unfinished 

charging demand at time k that must be completed by time t. 
With this, we define the state of system at time t as 

Dt = [It, dt, dj+i, ■ • ■ , dy], (5) 

where it is the base load at time t, d\, is the total unfinished 
charging demand at time t that must be completed by time 
f. Let represent the random arrival events at time t. is 
defined as 

it = [‘-t,VuvUw- ( 6 ) 

where it is the base load at time t, rf^, is the total charging 
demand that arrive at time t and must be fulfilled by time t', 
et is the latest deadline among the PEVs that arrive at time t. 
Then, the state transition, defined as 

Df+i := (7) 

is calculated as follows: 

It+l = H+l (8) 


and 




t' 

+- 

dlY = 

d\,- 

- E 



+4Y,'if = ,T. 


(9) 


Here, [a;]+ = max{x,0}. With the above definitions of^stem 
state and state transition, we are now ready to rewrite 0 into 
the following finite-horizon dynamic programming promem. 


Qk (Life) min f {sk -f L) [Q^+i idi^k^ Dfc, ^^+ 1 ))] 

(10a) 

T 

s. t. < Sk <'^ dt, (10b) 

t=k 


where Qk(Dk) is the optimal value of the dynamic program¬ 
ming at time k. The left side of (lObi ensures all charging 
demands to be satisfied before their deadlines. The right side 
of ( |10b[ ) implies that the total charging power up to a certain 
time cannot exceed the total demands that have arrived up to 
that time. By slight abuse of notation, in the rest of the paper 
we denote the optimal solutions to both the online and offline 
problems as s^, when no confusion arises. The actual meaning 
of will be clear from the context. Suppose that is the 
optimal solution to ( [T0| ) at stage k. Then, the total cost at the 
end of system time, denoted by T' 2 , is provided by 


T 

k=l 


( 11 ) 









Note that (lOai comprises nested expectations with respect 
to the random PEV arrivals in the future time stages. Except 
for few special cases, it is hard to provide the closed-form 
of the optimal solution to ([T§. On the other hand, 
can be solved by the standard numerical methods, such as 
backward reduction and the sample average approximation 
(SAA) based on Monte Carlo sampling techniques 
p2) . These algorithms typically incur a computational com¬ 
plexity that grows exponentially with both the time span T 
and the dimensions of state and decision spaces. Note that ( [T 0 | ) 
involves continuous state and decision spaces. Discretization 
of these spaces leads to a curse of dimensionality, rendering 
the computational complexity prohibitively high. 


III. MPC-based Online Charging Algorithm 


In view of the extremely high complexity of standard 
numerical methods, we are motivated to obtain a near-optimal 
solution by solving a much simpler problem, which replace 
all exogenous random variables by their expected values. This 
is referred to as the expected value problem 0 -® or the 
MPC approach p^-p3) in the literature. Notice that the first- 
order moment, i.e., expectation, of a random process is much 
easier to estimate than the other statistics, e.g., variance or the 
probability distribution. Thus, the assumption of knowing the 
expected values is weaker than the assumptions in other EV- 
charging algorithms | | 20 ) , which assume that the probability 
distributions of the random processes are known. 

Instead of solving the problem using generic convex op¬ 
timization approaches, we propose a low-complexity online 
Expected Load Elattening (ELE) algorithm by exploring the 
load flattening feature of the optimal solution to the expected 
value problem, as shown in Section III-A Section |III-B pro¬ 
vides the theoretical analysis of the performance gap between 
the optimal solution to the expected value problem and the 


optimal solution to (lOi. 


A. Algorithm Description 

Denote the expectation of as ,yLty], 

where vt = = E[77j,],Vf' = t, - ■ ■ ,T. Replacing 

in ( [T0| with /xt, we obtain the following deterministic problem: 

T 

min f{sk + lk)+ y'' f{st + i't) (12a) 

t = k + l 

3 3 3 3 

s. t. '^St>'^dt+ yy /x™,Vj = A:,-• • ,r, (12b) 

t = k t = k m=k + l n = m 

yyst<yyd*+ yy (i 2 c) 

t = k t = k m=k + l n = m 

In each time fc, we solve problem ( [T^ and obtain the optimal 
charging solution s^. Then, problem |l2[i is resolved with the 
updated according to the realization of the PEVs arrived in 
next time. So on and so forth, we obtain the optimal charging 
solution s*f. for time stage fc = 2, • • • , T. The total cost at the 


end of system time, denoted by T' 3 , is defined as 


> 1^3 =^/(4 + 4 ), ( 13 ) 

k = l 

where is the optimal solution to ([T 2 ji at time stage k. The 
solution to ( [T 2 | ) is always feasible to ^ 10 [) in the sense that 
it always guarantees fulfilling the charging demand of the 
current parking PEVs before their departures. This is because 
the constraints of Sk in ( fTO] ) are included in •Ill- 

Due to the convexity of /(•), the optimal solution is the one 
that flattens the total load as much as possible. By exploiting 
the load flattening feature of the solution, we present in 
Algorithm [T| the online ELE algorithm that solves ( [T^ with 
complexity O(T^). The online ELE algorithm have a lower 
computational complexity than generic convex optimization 
algorithms, such as the interior point method, which has a 
complexity 0 (T^-®) pT) . Notice that similar algorithms have 
been proposed in the literature of speed scaling problems pT) , 
p4| and PEV charging problems 0. The optimality and the 
complexity of the algorithm have been proved therein, and 
hence omitted here. The algorithm presented here, however, 
paves the way for further complexity reduction to 0 ( 1 ) in 
Section |IV] Eor notation brevity, we denote in the online ELE 
algorithm 

for f" = fc, • ■ • ,r,t'= fc, 

*" foit" = ,T,t' = k + l,--- ,T. ' 


The key idea of online ELE algorithm is to balance the 
charging load among all time slots fc, • • • , T. Specifically, step 
3 - 5 is to search the time interval [x* ,j*] that has the maximum 
load density among current time slots and record the maximum 
load density. The time slots with maximum load density are 
then deleted, and the process is repeated until the current time 
slot k belongs to the maximum-density interval, i.e., i* = k. 


Algorithm 1: Online ELE Algorithm 


input : T)k,Pt,t = fe + 1, • • • ,T 

output: Sk 

1 initialization i = 0 ,j = 0 ; 

2 repeat 

3 For all time slot i = k, - ■ ■ ,T, j = i, ■ ■ ■ ,T, compute 


* ,J 
Set 


arg max { + I't') ^ 

k<i<j<T^ J — 1+1 ’ 


y = 


j* — i* -I- 1 


(16) 


5 Delete time slot i*,■ ■ ■ ,j* and relabel the existing time 
slot t > j* as t — j* + i* — 1 . 

6 until i* = k\ 

7 Set Sk =y* - Ik- 


B. Optimality Analysis 

In this subsection, we analyze the optimality of the solu¬ 
tion to (12i. Notice that MPC approximates the non-causal 
















random variables by their expected values regardless of their 
distribution functions. As a result, such approximation may 
lead to unacceptably large performance loss, depending on the 
distribution of the random variables. Therefore, the MFC ap¬ 
proximation is not always justifiable. A well-accepted metric. 
Value of the Stochastic Solution (VSS) is adopted to evaluate 
optimality gap between the optimal online solution and the 
solution to the expected value problem Previous 

work, e.g., 0, |T§, mainly evaluates VSS using numerical 
simulations. Whereas in our analysis, we show that VSS is 
always bounded regardless of the distribution of the future 
EV charging demands. This provides a strong theoretical 
justification of adopting MPC approximation to solve our 
problem. 

Let S denote a scenario, which is defined as a possible 
realization of the sequence of random load demand p 2 ) , 

“ = [^2i ■ ■ ■ > ^t]- (1”^) 

Here, we treat as deterministic information since the 
demand of PEVs arrived at the first stage is known by the 
scheduler. Let $i, $2 and $3 be the expectation of the optimal 
value of the offline problem (|4i, the online problem ( [T 0 | ) 
and the expected value problem (^, respectively, where the 
expectation is taken over the random scenarios. That is, 

$1 = Es [®i(H)] , $2 = Eh [^ 2 (H)] , $3 = Eh [^ 3(H)] . (18) 

It has been proved previously 0 ^ 0 that 

'I’l < >f>2 < 'I>3. (19) 


To assess the benefit of knowing and using the distributions 
of the future outcomes, the VSS is defined as 

VSS = $3 - $2. (20) 

To show that the online ELE algorithm yields a bounded VSS, 
we need to bound $3 and $ 2 . Generally, it is hard to calculate 
$2 or analyze the lower bound of $2 directly 0, p0| . Thus, 
we choose to analyze the lower bound of $1 instead, since 
( [T^ shows that the lower bound of $1 is also the bound of 
$ 2 - In what follows, we will show the lower bound of $1 in 
Proposition and the upper bound of < 1>3 in Proposition 
Proposition 1: 

E:ii d\ + EL2 i^i + ELi 

T 

Proof: Please see the detailed proof in Appendix 

Let 0{t) be the set that 0{f) = {{m,n)\em > t,m = 
1, • • • ,t,n = t, ■ ■ ■ , Cm}- Then, we show that $3 is bounded 
by Proposition]^ 

Proposition 2: For any distribution of = I,-- - ,r, 
there is 


$3 ^ E 

1 _ 

E + A 


t=i \ 

^{m,ri)GO(t) J 


Proof: Please see the detailed proof in Appendix 

Now, we are ready to present Theorem 0 which states that 



$1 > Tf 


the VSS is bounded for any distribution of random variables. 

Theorem 3: For any distribution of random vector f = 
1, • • • ,T,n = t, ■ ■ ■ ,T, there is 


VSS < E 

1 _ 

E +1' A 


t=l \ 




( 23 ) 


where E = Etii 4 + EL 2 'LT=t + ELi 

Theorem 0 can be easily derived by Proposition 0 and 
Proposition 0 In practice, the performance gap between the 
online ELF algorithm and the optimal online algorithm is often 
much smaller than the bound of VSS. This will be elaborated 
in the numerical results in Section lYl 


IV. Online ELF Algorithm under First-order 
Periodic Random Processes 

Notice that the complexity of 0(T^} of online ELF algo¬ 
rithm mainly comes from step 3, which exhaustively searches 
the maximum-density period [i*,j*] over all subintervals 
within [k,T]. When the random arrival process is first-order 
periodic stochastically Q we argue that the searching scope 
can be reduced to one period from the whole system time 
T. Thus, the complexity of step 3 is limited by the length 
of a period instead of T. As a result, the complexity of 
the algorithm reduces from 0{T^) to 0 ( 1 ), implying that 
it does not increase with the system time T, and thus the 
algorithm is perfectly scalable. In practice, the arrival process 
of the charging demands are usually periodic stochastically. 
For example, the arrival of charging demands at a particu¬ 
lar location is statistically identical at the same time every 
day during weekdays (or weekends), the National Household 
Travel Survey (NHTS) 2009 gathers information about daily 
travel patterns of different types of households in 2009, and 
shows that the daily travel statistics (e.g.. Average Vehicle 
Trip Length, Average Time Spent Driving, Person Trips, 
Person Miles of Travel) are very similar for each weekday or 


weekend, but different between weekday and weekend 1281. In 


Section IV-A we investigate the case when the random load 
demand process is first-order periodic. In Section IV-Bl we 


provide a closed-form solution to (|I 2 | for a special case when 
the load demand process is the first-order stationary. 


A. First-Order Periodic Process 

In this subsection, we consider the case when the arrival 
process is first-order periodic. Specifically, the first-order pe¬ 
riodic process means that the first-order moment (i.e., mean) 
of the random process is periodic. That is, at current time stage 
k, for all f = /c -b 1 , • • • , T, we have 

Ft = E[lt] = Mt+p, (24) 

where is the random arrival events at time t, is the 
expectation of and p is the length of period. Then, instead 


*The first-order periodic stochastic process is defined as a stochastic process 
whose first-order moment, i.e., mean is periodic. That is, the mean of the 
random arrival events, i.e., is periodic. However, the actual realizations of 
the an'ival events are uncertain and not periodic. 











of considering for f = k + ■ ■ ■ ,T, we only need to 

consider for one period, i.e., for t = k + l,k + p: 

_ r fc + l fc + l fc + 1 pv pjI 

Mfe + 1 [^fc + 1 ? f^k + 1 ’ Mfe + 2 ?'■' 5 /^fc + ei’^’’'’ 


(25) 


_ r K-\-p Ac+p 

f^k+p — [I'k+p, Pis^p, P-k+p +1 


k-\-p 

>/^fe + e. 


, 0 ], 


Here, e„ < T, n = 1, • • • ,p is the maximum parking time 
for PEVs arriving at time k + n. Specially, we define e as 
e = max{e/c+i, efc+ 2 , • • • ,ek+p}- We decompose the search 
region {i,j\i = k,--- ,T,j = ),••• ,T} into three sub- 
regions, defined as Hi = {i,j\i = k, j = k, ■ ■ ■ ,fc-|-e}, 02 = 
{i,j\i = k,j = k + e + I, - ■ ■ ,T} and IIs = {i,j\i = k + 
1, • • • ,T,j = i, • • • , T}, respectively. We denote by X, Y, Z 
the maximum densities of region 01 , 112 , 113 , respectively. 
Indeed, the largest of X, Y, Z is the maximum density of the 
interval C [k,T] over all possible pairs i,j G {i = 

k,--- ,T,j = * *,••• ,T}. Let [ii, ji], [* 2 ,^ 2 ], [* 3 ,j 3 ] be the 
intervals with the maximum density over region 01 , 02 , 03 , 
respectively. By definition, zi = *2 = k. Similar to the 
stationary case, X can be calculated by searching ji over 
{fc, • • • , fc -I- e}. That is. 


X = max + Vn) + k‘m 

k<t<k + e n — fc + 1 

Moreover, Lemma shows that Y and Z can be calculated 
once the maximum density of interval [/c + 1 , fc + e] has been 
obtained. Lirst, we introduce some definitions which help to 
show Lemma[T] Let O4 be a subset of O3, where O4 is defines 
as O4 = {i,j\i = /c + 1, • • • ,k + e,j = i, - ■ ■ ,k + e}, and 
[i,j] be the interval with maximum density of region O4, i.e.. 






= arg max 

k-\-l<i<j<k-\-e 

Lemma 1: The maximum densities of O 2 and O 3 are cal¬ 
culated by 




j 2 -k+l is - *3 -I- 1 

respectively, where 13 = i, 

I max{i, fc-I-e-I- 1 }, if j<i+p, 

t ?2 ^ 1 ~ 

I i -I- (r — l)p, otherwise. 


and 


J3 


j, ifi<*+P, 

j + {r — l)p, otherwise. 


(27) 


(28) 


(29) 


Proof: Please see the detailed proof in Appendix 

Based on Lemma [T] we can modified the searching region 
of step 3 of online ELL algorithm as follows: 

• if j < i + p, the interval with the maximum density 

during time stages [fc -f 1,T] is Then, in step 3 
of the online ELL algorithm, the search region of i,j 
is reduced from {z,i|z = k,--- ,T,j = z, ••• ,T} to 
{z,_j|z = fc, ••• = ,ij}. 

• If i > i+P, the interval with the maximum density during 
time stages [k -f 1,T] is [i,j + {r — l)p]. Then, in step 


3 of the online ELL algorithm, the search region of i,j 
can be reduced from {i,j\i = k, - ■ ■ ,T,j = i, - ■ ■ ,T} 
to {ij\i = k,--- ,i,j = ... + (r- l)p}. 

As a result, the searching region of the online ELL algorithm 
is only related to [fc -f 1, /c -f e] instead of T. Thus, the 
computational complexity of the online ELL algorithm is 0(1) 
instead of O(T^) under first-order periodic process. 


B. First-order Stationary Process 


In this subsection, we show that the optimal solution to 


( 121 can be calculated in closed form if the arrival process is 
first-order stationary. Due to the page limit, we only provide 
the main results here. By first-order stationary, we mean that 
the statistical mean of i.e., vt and /r*,,f' = f, • • • ,T only 
depends on the relative time difference t = t' — t, but not the 
absolute value of t. We can then replace vt by v and replace 
p\, by Pp, where t = t' — t. Then, is no longer a function 
of t, and can be represented as 

P= ,/Xe,0, ,0], (30) 


where e is the maximum parking time of a PEV. We denote 
by X,Y,Z the maximum densities of region ni,n 2 ,n 3 , 
respectively. Then, X is calculated by 


X = 


max 

A3<n<fc+e 


EEfc + E)=i i'n — k — j l)pj Ik — I 

n — k 1 


■ -f v 


(31) 


and y, Z are provided in Lemma 

Lemma 2: The maximum densities of n 2 and 03 are 
achieved by setting 12 = k,j 2 = T, 13 = fc -I- l,j 3 = T, 
and calculated by 


E"i" + E"if (T -k-j + l)p, +lk-u 
- T^kTi -+ 

j:-t!iT-k-j + i)p, 

-* J_ 1 / 


(32a) 

(32b) 


Proof: Please see the detailed proof in Appendix The 
largest of X, Y, and Z is the maximum density of the 
interval [z*,E] ^ [k,T] over all possible pairs i,j G {i = 
k, - ■ ■ ,T,j = z, • • • , T}. Specially, if X or L is the largest 
one, then k is already contained in the maximum-density 
interval, and thus X or L is the optimal charging rate at time 
k. On the other hand, if Z is the largest, then the maximum- 
density interval, i.e., [fc-|-l,T], does not include k. Lollowing 
Algorithm [T] we will delete the maximum-density interval and 
repeat the process. Now, time slot k is the only remaining time 
slot after deletion. This implies that all charging demands that 
have arrived by time slot k should be fulfilled during time slot 
k. These arguments are summarized in Proposition which 
provides the closed form solution to (HI- 

Proposition 3: When the random load demand process is 
first-order stationary, the optimal charging schedule to ([T2|i is 











given by the following close-form: 


x-h, 

if X = max{A:, V, Z}, 

(33a) 

Y-h, 

if Y = max{X, y, Z}, 

(33b) 

k-\-e 



t^k 

otherwise. 

(33c) 

V. 

Simulations 



In this section, we investigate the performance of the 
proposed online ELF algorithm through numerical simulations. 
All the computations are solved in MATLAB on a computer 
with an Intel Core 13-2120 3.30GiJz CPU and 8 GB of 
memory. For comparison, we also plot the optimal solution to 
which is obtained by SAA method | |22| , and a heuristic 
solution by the online AVG algorithm Q, which is obtained by 
charging each PEV at a fixed rate, i.e., its charging demand 
divided by its parking time. Let the expected cost of AVG 
algorithm is denoted by $ 4 . Define the relative performance 
loss of ELF and AVG compared with the optimal online 
solution as and respectively. Similar to 0 

Q, we adopt an increasing quadratic cost function in me 
simulations, i.e., f{st + h) = {st + hY■ Note that the cost 
function is increasing and strictly convex, since the load St + lt 
is always non-negative. 

A. Average Performance Evaluation 

In this subsection, we evaluate the average performance of 
the online ELF algorithm under three different traffic patterns, 
i.e., light, moderate, and heavy traffics. In particular, the 
system time is set to be 24 hours, and each time slot lasts 
10 minutes. The PEV arrivals follow a Poisson distribution 
and the parking time of each PEV follows an exponential 
distribution @- 1 ^. The mean arrival and parking durations 
of the three traffic patterns are listed in Table The main 
difference lies in the arrival rates at the two peak hours, i.e. 
12 : 00 to 14 : 00 and 18 : 00 to 20 : 00. The settings of the 
peak hour match with the realistic vehicle trips in National 
Household Travel Survey (NHTS) 2009 p 8 ) . Specially, the 
average number of total PE Vs simulated in scenario 1, 2 
and 3 are 104,204 and 304, respectively. We choose the 
base load profile of one day in the service area of South 
California Edison from |7). Each PEV’s charging demand is 
uniformly chosen from [25, 35]A:lU/i. Each point in Fig. 
Fig. § and Table [n| is an average of 10® independent instances 
of the scenarios listed in Table In Fig. [^ the total loads 
s*f. are plotted over time. We notice that the total load of 
the online ELF algorithm follows closely to that of optimal 
online algorithm, whereas that of the AVG algorithm has a 
larger gap from the optimal online algorithm. The average 
costs normalized by that of the optimal online algorithm 
are plotted in Fig. Moreover, the VSS and the relative 
performance loss are listed in Table [I^ Both the figure and 
the table show that ELF performs very close to the optimal 
online algorithm. The VSS and the relative performance loss 
are no more than 0.1536 and 0.38% respectively. In contrast, 
the relative performance loss of AVG algorithm is up to 5.82%, 


TABLE I 

Parameter settings oe the PEV traffic patterns 


Time of Day 


Arrival Rate (PEVs/hour) Mean Parking 


S. 1 S. 2 


S. 3 


Time (hour) 


08 : 

: 00- 

- 10 : 

: 00 

7 

7 

7 

10 

10 : 

: 00- 

- 12 : 

: 00 

5 

5 

5 

1/2 

12 : 

: 00- 

- 14 : 

: 00 

10 

35 

60 

2 

14 : 

: 00- 

- 18 : 

: 00 

5 

5 

5 

1/2 

18 : 

: 00- 

- 20 : 

: 00 

10 

35 

60 

2 

20 : 

: 00- 

- 24 : 

: 00 

5 

5 

5 

10 

24 : 

: 00- 

- 08 : 

: 00 

0 

0 

0 

0 


TABLE II 

Average performance comparison under three trafeic patterns 


Scenario 

VSS 

^3 — *^*2 

^>2 

' i '4 —‘ i >2 

^2 

1 

0.1178 

0.19% 

3.50% 

2 

0.1319 

0.28% 

4.46% 

3 

0.1536 

0.38% 

5.82% 


which is more than 15 times of that of ELF algorithm. The 
relative performance loss of the approximate online algorithm 
reflects the percentage of the extra cost compared with the 
optimal online algorithm. Obviously, the performance loss is 
always the smaller the better. For example, from the report 
of Rocky Mountain Institute, the average electricity rate is 
11.2cents/kWh, and the average load for a charging station is 
lOOkW p^ . Then, the expected electricity cost of a charging 
station for one years is $967680. 6 % relative performance loss 
means AVG algorithm leads to $58060 extra cost, while the 
proposed online algorithm with 0.38% relative performance 
loss leads to $3677 extra cost, or an annual saving of $54383 
compared with the AVG algorithm. 

B. Complexity of The Online Algorithm ELF 

In this subsection, we verify the computational complexity 
of the online ELF algorithm and also compare the complexity 
of online ELF algorithm with that of the optimal online 
algorithm and online AVG algorithm, respectively. We still 
adopt the SAA method as the optimal online algorithm. Since 
the complexity of SAA is too high, truncation is often adopted 
in the SAA to reduce the complexity at a cost of performance 
loss Q. As such, we also simulate the complexity of truncated 
SAA with a truncation period of 3 hours. Each PEV’s charging 
demand is uniformly chosen from [25,35]fclU/i, the arrivals 
follow a Poisson distribution and the parking time of each 



Fig. 1. 


Normalized costs of three algorithms in three scenarios. 
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(a) Scenario 1: light traffic 


(b) Scenario 2: moderate traffic 


(c) Scenario 3: heavy traffic 


Fig. 2. Base load and total load of three algorithms in three scenarios. 



Fig. 3. CPU computational time over the system time. 

PEV follows an exponential distribution, where the arrival rate 
and mean parking durations are the same as that in the peak 
hours of scenario 2 in Section V.A. We simulate 10 cases in 
total, where the system time are set to be 1, 2, • • • , 10 hours. 
For each algorithm, we record the CPU computational time at 
each time stage and calculate the average CPU computational 
time as the sum of CPU computational times of all time 
stages divided by the number of time stages. Each point in 
Fig. is an average of 100 independent instances. Fig. 
shows that the CPU computational time of both the online 
algorithm ELF and AVG almost grow linearly with the system 
time. In contrast, for the SAA based online algorithm with 
or without truncation, the average CPU computational times 
grow very quickly as system time increases. We notice that 
when system time increases to 4 hours, the optimal online 
algorithm without truncation consumes more than 2 hours and 
the optimal online algorithm with truncation consumes more 
than 30 minutes. Meanwhile the proposed online algorithm 
ELF only takes several seconds. It is foreseeable that the 
computational complexity of optimal online algorithm will 
become extremely expensive as we further increase the system 
time. 

C. Performance Comparison with Online Algorithm OR¬ 
CHARD 

In this section, we compare the proposed online algorithm 
ELF with the online algorithm ORCHARD proposed in refer¬ 



Fig. 4. Load comparison of ORCHARD and ELF 

ence |j^ on the properties of the optimality and computational 
complexity. First, we evaluate the average performance of the 
proposed online algorithm and online algorithm ORCHARD. 
To facilitate the comparison, we also adopt the average per¬ 
formance of optimal online algorithm as a bench mark. For 
the system parameters, we use the default settings of scenario 
1 in Section V.A. We simulate 10® cases and plot the total 
load (the sum of the base load and the PEV charging load) 
over time for the three online algorithms in Fig. In addition, 
the average performance ratios normalized against the optimal 
online solution are shown in Fig. [^respectively. Fig. [^ shows 
that compared with online algorithm ELF, the online algorithm 
ORCHARD always produces a larger gap from the optimal 
online solution. From Fig. [^ we can see that the proposed 
online algorithm ELF achieves a much lower expected cost 
than the online algorithm ORCHARD, which indicates that 
online algorithm ELF owns a better average performance than 
online algorithm ORCHARD. 

To compare the computation complexity of online algorithm 
ORCHARD and ELF, we adopt the similar case study in 
Section VD of reference Q. Specifically, we simulate the 
CPU computational time of online algorithms by varying the 
arrival rates of the PEVs during one day. For the system 
parameter, we use the same settings as scenario 1 in Section 
V.A except the arrival rates, which are assumed to be the same 
during 8 : 00 — 18 : 00 and 0 after 18 : 00. We vary the arrival 








































Fig. 5. Average performance ratios of online algorithm ORCHARD and 
ELF 



Fig. 6. CPU computational time over the number of PEVs. 

rate in 8 : 00 —18 : 00 from 10 to 50 (PEVs/hour) that leads to 
the average number of total PEVs during one day varies from 
100 to 500. Eor each specified average number of PEVs, we 
simulate the average performance of 10 ® independent instances 
for the online ELE algorithm, the optimal online algorithm and 
the online ORCHARD algorithm, respectively, and record the 
average CPU computational times of each cases for the three 
algorithms, respectively. The results are plotted in Pig.|^ Fig |^ 
shows that the average CPU computational time of the optimal 
online algorithm grows quickly as the number of total PEVs 
increases, while the average CPU computational time of the 
online ELE algorithm and the online ORCHARD algorithm 
grows slowly as the number of total PEVs increases. When 
the number of PEVs is 200, the average CPU computational 
time of the optimal online algorithm without truncation is 
more than 24 hours. Even for the the optimal online algorithm 
with truncation, the average CPU computational time is about 
100 minutes. Whereas, the proposed online algorithm ELE 
only takes about 4 minutes. Eig. also indicates that the 
computational time of online ELE and online ORCHARD as 
the number of the PEVs increases. 

As a conclusion, the case study shows that the online 
algorithm ELE has a better average performance than online 
algorithm ORCHARD. It also indicates that the CPU compu¬ 
tational time of online ELE and online ORCHARD are similar. 

VI. Conclusions 

In this paper, we formulate the optimal PEV charging 
scheduling problem as a finite-horizon dynamic programming 


problem. Instead of adopting the standard numerical meth¬ 
ods with high complexity, we provide a MPC-based online 
algorithm with 0(T®)-complexity. We rigorously analyze the 
performance gap between the solution of the MPC-based 
approach and the optimal solution for any distribution of 
exogenous random variables. Moreover, we show that the 
proposed algorithm can be made scalable under the first-order 
periodic process of load demand. Besides, our analyses are 
validated through extensive simulations. 

Appendix 

A. Proof of Theorem 

We use the inductive method to show that through EDE 
scheduling, all the PEVs can be fulfilled charging before 
deadlines. Eor n = 1, 0 becomes 

di> Si> di. (34) 

Thus, by EDE scheduling, we can first satisfy the demand of 
PEVs whose deadline at time stage 1. That is, for any PEV 
i G = 1}, we set 

^ii — di. (35) 

Assuming that for all time stage m, EDE scheduling can fulfill 
charge all the PEVs which depart at or before time stage m, 
i.e., there exists at least a set of Xifs that satisfy 

df'’ 

^ Xit = di^i G < m}, (36a) 

t=&^ 

Xit > 0, Vf = G < m}. (36b) 

Since 

m m 

E (37) 

then, Y.7=i St - J2T=i di represents the amount 

of power which is outputted from the charging station during 
time stage 1, • • • ,m and charged to the PEVs with deadline 
after time stage m. By EDE scheduling, once the PEVs which 
depart at time m have been fulfilled charging, we will first 
charge the PEVs which depart at time stage m+ 1. Thus, if 

m m 

E®*-E E E d„ (38) 

=^t} 

we finish charging of PEVs with deadline m + 1, and then go 
to charge the PEVs with deadline m + 2. If 

m m 

E®*-E E E d^, (39) 

then the PEVs with deadline m + 1 have been charged as 
power YT=i St-J2T=i m + 1. 










Since 

then, 

Sm+l ^ 


m + 1 m + l 

t ^ I] 

‘=1 ie{i|t?'>=t} 


(40) 


/ 


i6{i|t 


E 


di — 


.+ 1 } 


E - E E ’ (41) 

‘=1 ie{i|t('>=t} 


which means all the PEVs with deadline m+l can be fulhlled 
charging. This is because we will charge the PEVs with 
deadline m + l first by the EDE scheduling. Thus, there exists 
at least a set of Xu's that satisfy 


1) If x*f^ = 0 for a particular PEV i at a time slot ti G 
{fthen, by complementary slackness, we 
have ujiti > 0. Erom ( [47a| , 

/^(sji + (ji) = Ai + Witi. (48) 

2) If x*f^ > 0 for PEV i during a time slot t 2 G 

we can infer from (47ci that ujn^ = 0. 

Then, 

f(st,+lt,) = X,. (49) 

On the other hand, since f{st+lt) is a strictly convex function 
of St +lt, then f'{st + lt) is an increasing function. Erom the 
above discussions, we get the following two conclusions: 




E ^ = m + 1}, 


Xt,m +1 > 0,Vi e = m + 1}. 

Combining ( [36l l and ( |42| ), we get that all the PEVs whose 
deadline at or before stage m + l can be fulhll charging, i.e., 
there exist at least a set of Xu’s that satisfy 


(42a) 

(42b) 


1) If x*^^ > 0,x*^^ > 0, then by ( p9] l, 

f'isti + hi) = f'{st^ + ^* 2 ) = Ai. (50) 
Due to the monotonicity of we have + hi = 

5(2 + ((2 ■ 

2) If xh^ = > 0, then by (48i and (491, there is 


4=) 


E = di,\/i G < m + 1}, 


(43a) 


As) 


Xit > 0,Vf = G < m + 1}. (43b) 

Therefore, we can conclude that by EDE scheduling, there 
always exists at least a set of Xu’s that is feasible to (|^. This 
completes the proof. ■ 

B. Proof of Theorem ^ 

Eirst, we show that if there exists a PEV parking in the 
station at both time ti and t 2 , i.e.. 


f'{sti + hi) — Ai + ujit^ > f'{st2 + (( 2 ) — Ai. (51) 

Since f'{st) is a increasing function, we have +hi > 
5(2 + h^ ■ 


Consider two function f{st + h) and f{st + h)- Let x*f. and 
xh denote the optimal solutions to (j^ with f{st + h) replaced 
by /(5( + h) and f{st + h), respectively. Dehne Sj, Sj as 


= E E = h 

iGT(t) iGT{t) 


T 


(52) 


respectively. Suppose that there exists a time slot ti such that 


5,. < 5, 


(53) 


ti,t2 G lE, • • • ,E}, 

(44) 

Since 

T T 


and 

x*ti > 0, x*t2 > 0, 

then the optimal total loads at time ti and ^2 

5(i T hi ^ 5(2 + ((2 • 

(45) 

must satisfy that 

(46) 

(=1 (=1 

there must exist another time slot t 2 such that 

-:ic 

St2 > St2 

(54) 

(55) 

The Karush-Kuhn-Tucker (KKT) conditions to the convex 
problem (|^ are 

/'( E Xit + k) - Xi - Xit = 0,i G J\f,t = 

and 

Thus, we can 

1 —* 1 —* 

5(1 + 5(2 = 5t^ + S(2 

find a PEV i G Af such that 

(56) 

i£I(t) 

(47a) 


Xiti < Xt^,X,t^ > x,t^. 

(57) 

(<=) 

Xi{di— E Xit)^0,iGN', 

(47b) 

As a result. 

* 

to 

V 

0 

(58) 

XitXit = Q,i G N ,t = \ ■ ■ ■ , , 

(47c) 

since x^^ > 0 

. Based on (46 1, there is 

5(2 + ((2 < 5(^ + hi ■ 

(59) 


where A, w are the non-negative optimal Lagrangian multipli¬ 
ers corresponding to ( [2b] i and ( [2c| ), respectively. We separate 
our analysis into the following two cases: 


Combining ([53)([56)([59]l, we get 

5(2 + ^(2 < 5(2 + ((2 < S(^ + hi < Ki + hi ■ 


(60) 













Since f{st + h) is a strictly convex function of St + It, then, 
based on and ( |60l l, we have 

f {Ki + 1 ‘ti) + f {st^ + k 2 ) > f {sti + ki) + f {st^ + lt 2 ) ■ (61) 


J 2 t= 2 '^T=t ^^t■ Theiij by Jensen’s inequality |30|, 




I't 


T 


This contradicts with the fact that the Sj is the optimal 
total charging rate for objective function f{st + It)- There¬ 
fore, the optimal charging solution Sj is the same for any 
strictly convex function f{st + U)- Next, we show that opti¬ 
mal solution St is a load flattening solution that minimizes 

E^i(st + It — subjecting to ( |4b| ) and ( |4cl i. 

Based on the argument that Sj is the same for any strictly 
convex function f{st -f h), then it is equivalent to show that 
J 2 t=ii^t + h — is a strictly convex function of 

St + It- Since 




T 

which indicates that 


T 


(62) 




is a constant. Then, we see 


^ /(st(H"') -I- It) 


t=i 

T 


(63) 


+ lt) -I- (1 — A)^/(st(H") + lt) 


for all A S [0,1]. On the other hand, let be the optimal 

solution that minimizes E^i fi^t + h) subject to (4bi - (4c I 
with H = s'”. Then, 


^/(4(h'") + /0<E/(«*(“") +^‘)- 


(64) 


Combining ([63) and (64i, we have 


t=i 

T T 

<aE /(®t(“0 + ^t) + (1 “ '^fi^ti'^") + ^t)- 


(65) 


Thus, we have established the convexity of over the 

set of H. Therefore, we have 

E 1^1(3)] > ^i(E[E;]). (66 ) 

On the other hand, based on the definition of St,dl,iil,i = 
!,■■■ ,ei,j = t,--- ,et, we have ELi + 


= Tf 


EEirfl + E? 




^t 


T 


(67a) 


(67b) 


This completes the proof. ■ 

D. Proof of Proposition [^ 

St < holds for all stage t. On the other hand, we 

have 


n=t n=t — l n=t n=t 


( 68 ) 


that E^i('5* + h— ^ strictly convex function 

of St + h- This completes the proof. ■ 

C. Proof of Proposition [^ 

First, we show that 'I'i(H) is a convex function of H. For 
any H, define S((H) as the optimal solution that minimizes 
'I'i(H) subject to ( |4^ - (j^. For an arbitrary pair of 2' and 
2", let s'” = AH'-|-(1 —A)H"VA S [0,1]. Then, the charging 
schedule such that = Asj (H')-l-(l —A)S( (H") 

still satisfies (4bi - (4ci with H = S'” due to the linearity of 
the constraints. Based on the convexity of f{st + It), we have 


where the inequality holds since st_i > d\_\ and M. = 
{m|em, > t,m = 1,--- ,t} Thus, we have st < 

E(m.n)eO(t) ^ bounded set for t = 


Et=l fiJ2{m,n)GOit) Vn + ''*) 


1, • • • , T. Therefore, E 

upper bound of $3. Thfs completes the proof. 

E. Proof of Lemma 

We provide the proof by discussing the following two cases: 

1) If 3 > i + P, which means that \i,j\ and \i + p,j + p] 

overlaps with each other, then the density of interval \i,j -Pp] 
is higher than that of [i,j], and the density of [i,j -P 2p] is 
higher than that of [i,j -fp]. So on and so forth. Finally, we 
see that the interval [i,j -P {r— l)p] has the maximum density 
over region 113. Thus, we have = i,j 3 = j + {r — l)p, and 

'^j + (r-l)p(spk+er, n , ^ 

^ y^m—n > n) ( 69 ) 

j - 1 - (r - l)p - i -I- 1 

Likewise, for the region {i,j\i = k,j = fc -|- e -f 1, • • • ,T}, 
we have j 2 = j + {r — l)p, and 

,,n , \ 

T/- _ Z-^n=k _ \Z-^m=n ^ 70 ^ 

j + (r - l)p - fc + 1 • ^ ^ 

2) If j < i + p, then, the density of interval [i,j] is higher 
than that of [i,j -f p], and the density of [i,j -f p] is higher 
than that of [i,j -|- 2p]. So on and so forth. Finally, we see 
that the interval [i,j] has the maximum density over region 
113. Thus, we have 13 = i,j 3 = j, and the corresponding 
maximum density 

EE(E:E;M™+^n) 


z = 

i - * +1 

For the region 112, if j < k + e + 1, then 




S~-k + e+lfspk + e„ n , x 
Z^n = k \Z-^m=n > n) 


•\k-\-en 

jm=rL 

e-f 2 


If J > fc -f e -|- 1, then 

EL.(E:;EM- + ^n) 


Y = 


j - fc -I- 1 


(71) 


(72) 


(73) 


This completes the proof. 























F. Proof of Lemma ^ 


we have 


P{fj) = 


Let p{i,j) denote the maximum density of For any 
i = k,j = k + e + 1, ■ ■ ■ ,T, the density of interval [i, j] is 
given by 

SCi ECf d1 + h + U - k)u 

j - k + l 

^ SCi jj-k + l- f)pt + d’l + Ik-V 

j-k + l 

(74) 

To prove that the maximum density is achieved by setting 
j = T, we only need to show p{i,j) is a non-decreasing 
function of j for each given i, i.e., 


p{i,j) < Pii,j + l),Vfc + e + 1 < j < T - 1. 

Since 

I2i=i U - k + l- f)pt + ^ ^ 

•' t=i 

k + e+l<j <T -1, 

we have 

p{fj + 1) 

E|ii (i - fc +1 - t)p, + d? + Pt + h-^ 

j — k + 1 

^ EEi" {j-k + l-t)p^ + E?Z" d? + h-a 


(75) 


(76) 


-f a 


j - k 


+ V 


=Pihj), 

(77) 

which implies (751. Hence, Y is the maximum density of 
[k,j],j = k + e + 1, ■ ■ ■ ,T. Next, we show that Z is the 
maximum density of [fc + 1, T]. For any fc + l<j<j<T, 
the density of interval [i, f\ is given by 

E i-i+l s^i-i+2-k 
k=l 2^t=l 


P{fj) = 


Pt 


+ V 


3-i + l 


(78) 


+ V. 


3- i + l 

To prove that the maximum density is achieved by setting 
i = fc + l,j = T, we only need to show p{i,j) is a non¬ 
decreasing function of j for each given i, i.e., 

p{i,j) < p{i,3 + l),Vfc-f 1 < z < i < T - 1, (79) 

and a non-increasing function of i for each given j, i.e., 

p(i,3) > P{i + l,j),Vfc -fl<z-fl<j<r. (80) 

On one hand, since 

— z + 2 — t)pt yh 4- ^ < ■ < ■ < T 

- < > MZ) vfc + 1 < z < J < i , 


j - z + 1 


( 81 ) 


J -1 + 2 

^ EEr( 7 -z + 2-^)M. ^ 

3 - i + l 

= p{i,j), 

which implies ( |79l l. On the other hand, as 

Y.l=i{3 - i + ^ - -t)Pt 

7 - * 

then 

p{i + l,3) = 


(82) 


< ^ fit, Vfc-I-1 < z < j < T, 

(83) 


ELi(7-* + 1-0^4 




j - Z + 1 

= P{i,3), 

which implies ( |80l l. This completes the proof. 
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